9 dependence of CP 9 model 
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We apply to the CP 9 model two recently proposed numerical techniques for simulation of systems 
with a 6 term. The algorithms, successfully tested in the strong coupling limit, are applied to the 
weak coupling region. The results agree and errors have been evaluated and are at % level. The 
results scale well with the renormalization group equation and show that, for CP 9 in presence of a 
(9-term, CP symmetry is spontaneously broken at 9 — n in the continuum limit. 
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Understanding the role of the parameter in QCD and 
its connection with the "strong CP problem" is one of the 
major challenges for high energy theorists 0. Unfortu- 
nately euclidean Lattice Gauge Theory, our main non 
perturbative tool for QCD studies, has not been able to 
help us for two reasons: i) the difficulties in defining the 
topological charge operator on the lattice [2| and ii) the 
imaginary contribution to the action that comes from the 
0-term prevents the applicability of the importance sam- 
pling method. 

Concerning the 0-term effect on the low energy dy- 
namics, two dimensional CP^ -1 models are very inter- 
esting since they have #-vacua, a well defined topolog- 
ical charge and many features in common with four- 
dimensional Yang-Mills theories: they are asymptotically 
free, possess instanton solutions, and admit a 1/A expan- 
sion. Of course they share point ii) with QCD and this 
means that, from the point of view of numerical calcula- 
tions, the complex action problem is as severe as in QCD. 
Unfortunately the algorithms used in the past to tackle 
the complex action problem are limited to moderate lat- 
tice sizes and the fate of CP symmetry at = tt is still 
an open question 0,0. 

Recently we proposed two independent methods for 
simulating models with (9-term. Both are based on stan- 
dard (i.e. with real action) simulations of the model at 
imaginary 8 = —ih but then they use completely differ- 
ent strategies. The first is based on the reconstruction of 
the probability distribution function of the order parame- 
ter |5[ , the second introduces a particular transformation 
that allows to reconstruct the order parameter via an ex- 
trapolation of a suitably chosen quantity 0. Both can 
in principle suffer from systematics but the agreement 
of the two gives, a posteriori, a very stringent indication 
that systematic effects are indeed under control. 

We tested the two algorithms in some models where the 
analytical solution is available, and also a perfect agree- 
ment between the two approaches was found for CP 3 in 
the strong coupling region. The aim of present work is 
to extend this study to the continuum limit. 



Moving to the weak coupling region is a tricky issue. 
QpN-i moc [ e i s nave the problem of unphysical configura- 
tions characterized by topological structures of the size of 
one lattice spacing (called dislocations) 0, Q • The effect 
is relevant for small N while, for N — > oo, is suppressed 
and the model becomes a system of N free particles and 
N free antiparticles where 9- vacua are degenerate in en- 
ergy lid ] • If we are interested in the continuum limit 
the effect of dislocations has to be negligible not to spoil 
the scaling properties of the theory. For small N this 
forces to go to very large values of the correlation length 
and, consequently, untenable large volumes are manda- 
tory. 

Previous studies of CP^ -1 models at 9 = reported 
that TV = 4 is not large enough to suppress dislocations 
fll| . Indications of scaling of the topological susceptibil- 
ity Xt(9 = 0) were observed but violations of universality 
were also present: simulations with different actions lead 
to distinct asymptotic values 0, . It is not the scope 
of present work to analyze the contribution of disloca- 
tions to the topological charge at finite 9, and it is not 
easy to address this problem within our framework since 
different definition of the charge operator (like the one 
based on local polynomial definition or cooling) can not 
be directly applicable to 9 = —Hi simulations. What we 
know is that increasing N the situation is less severe and 
indeed, using the results of Campostrini et al 0], we 
can affirm that N = 10 is large enough. Therefore we 
concentrated our efforts on the scaling region of the CP 9 
model. 

We have adopted for the action the "auxiliary U(l) 
field" formulation 

S g = N/3 / ](z n+ll z n U nili + 

Un, M - 2) (1) 

where z n is a A-component complex scalar field that sat- 
isfies z n z n — 1 and U n , a is a [7(1) "gauge field". The 
topological charge operator is defined directly from the 
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17(1) field: 

p 

where U p is the product of the U(l) field around the 
plaquette and — ir < \og(U p ) < n. 

For both methods the action to simulate S = S g + 
Se = -ih is real and local. This implies that the main 
computational cost is practically equivalent to standard 
simulations of CP^ -1 model at 6 = 0. This is important 
since, as will be clear in the following, large volume sim- 
ulations were mandatory to have good quality raw data. 

We refer to the original papers for more information 
and summarize here the main points of the new methods. 

The first one [UEl reconstructs the p.d.f. of the topo- 
logical charge density that is introduced to define the fi- 
nite volume partition function as a discrete sum over the 
different topological sectors: 

Z v (6) = e- y l v ^e i9Vx . (3) 

X 

where x = n/V is the density of topological charge (n is 
an integer that labels sectors with different topological 
charge) and e~ v ^ v ^ is, up to a normalization constant, 
the above mentioned p.d.f. . From the measurement of 
the mean value of the density of topological charge as a 
function of h we use the saddle point approximation to 
get, up to 0(1/V) corrections, a numerical evaluation of 
fv(x), the derivative of fv{%) respect to x. After a fit of 
such derivative, an analytical integration is performed to 
reconstruct the desired probability distribution function. 
At the very end a multi-precision code is used to compute 
the sum (0 and, consequently, all the relevant physical 
quantities. 

The fitting procedure intrinsic in our approach avoids 
the pathological flattening behavior of the free energy as 
already noted by other groups 0. It is also important 
to notice that the multi-precision step needed to com- 
pute the observables (free energy density and topological 
charge density) from the definition © is not a bottleneck. 
As reported elsewhere 01 it can be shown that, at large 
volume, only the values of x such that fv{x) < 9v(^) 
(where gy (&) is the exact free energy density normalized 
to zero at 9 = 0) are relevant for the computation of the 
sum (J3J). This has some positive implications: 

(i) it is possible to fit only the relevant part of the data 
having very good fits with a suitably chosen function de- 
pending on few parameters; 

(ii) the number of terms relevant for the correct evalua- 
tion of the sum is much smaller than V and diminishes 
increasing j3 (recall that gv{6) goes to zero as (3 — > oo); 

(iii) the number of significant digits necessary for the 
calculation of J2J is related to the quantity e~ Vgv ^ (the 
smallest value the partition function takes) and, at fixed 
volume, the required precision in the multi-precision step 



diminishes when going deeper in the weak coupling re- 
gion. 

To take advantage of it we use an iterative procedure: 
from a first "trial fit" done with all the data available we 
can have an estimation of fv{%) and, using the multi- 
precision code, an estimation of <jy(7r). Using this ap- 
proximate information we can extract a trial value x max 
such that <jy(7r) < fv( x max)- Having a better deter- 
mination of the interval actually necessary for current 
analysis we can repeat the procedure steps. In practice 
the convergence is very fast (two, maximum three steps 
are enough in all cases) and, at the end, allows to fit only 
the relevant part of the data having very good fits. 

The second method |g uses x(h), the topological 
charge density, in a different way. It is possible to 
rewrite the sum (0 as an even polynomial in the vari- 
able z = cosh -| and to define the derivative of the free 
energy density respect to logz 2 as y(z) — x(h)/ tanh \. 
Using the transformation y\(z) — y(e*z) and plotting 
yx/y against y one gets typically a smooth function for 
small y 6] . In addition, in the weak coupling scaling re- 
gion, the function y\/y can be computed from numerical 
simulations at real h until values of y very close to 0. 
Hence we can expect a simple extrapolation to y — > 
(i.e. in the region corresponding to real 9) to be reliable 
and thus to obtain x(9) (of course the result has to be 
independent of the specific value of A). 

In the same spirit the effective exponent 7a = 
j^\og(yx/y) as a function of y will give the dominant 
power of y(z) as a function of z near z = or, equiva- 
lently, the behavior of x(ff) for 9 — > tv. If n j\(y — > 0) = 1 
CP symmetry is spontaneously broken at 9 — n, values 
between 1 and 2 indicate a second order phase transition 
with a divergent susceptibility and so on. Clearly the 
method is based on the feasibility of a simple extrapola- 
tion procedure. This can only be considered a working 
hypothesis but it has already been successfully tested in 
different models [fj. 

For both approaches systematic effects can be present: 
0(1/ V) corrections to saddle point solution and the ar- 
bitrariness of the fitting function in the first case (see 
the relative discussion in [5j), the reliability of extrap- 
olation in the second one. What is crucial is that the 
consistency of results coming from the two procedures is 
a strong indication that systematic effects are negligible 
and we are actually reconstructing the true 9 dependence 
of the model. 

For N = 10 and h = the scaling window starts at (3 
around 0.75 (see 0|) and we decided to concentrate our 
efforts from this point on. The first step is to establish 
the appropriate lattice sizes. From h = simulations 
at fixed /3 one realizes that the asymptotic value of the 
topological susceptibility is reached when the lattice side 
L is roughly larger than 8 times the correlation lengths. 
Using available data we get L > 40 for (3 = 0.8, 
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FIG. 1: The symbols are j x (y) in a L = 100, (3 = 0.80 lattice 
for A = 1.0. The continuous line is a quadratic fit to the data 
(reduced \ 2 is 1-20). 

L > 70 for (3 = 0.9 and so on. Convergence of \t to the 
infinite volume limit is even faster at imaginary 9 and 
consequently we expect the results at real 9 extracted 
from h ^ simulations will scale too. 

For h > simulations we considered different h val- 
ues (~ 50 where typically enough) in steps of 0.2. We 
are interested in precise measurements of the topological 
charge density as a function of h but, for each configu- 
ration, this quantity is a discrete number that takes the 
values n/V, where n is an integer between and V. It 
turns out that values of x much smaller than 1/V are 
very difficult to measure and, with affordable statistics, 
they suffer huge relative errors. The only way to have 
reliable measurements of x(h) for the smaller values of h 
is to consider large enough volumes. How large can be 
estimated using a linear extrapolation around h = 0: if 
we write x(h = 0.2) ~ \xt{h = 0) and impose the condi- 
tion x(h = 0.2) > 1/V we can get the minimum V value 
to be considered. Using again the available data [ll| we 
get that for (3 — 0.85, L has to be at least 100 while for 
(3 = 0.9 aL~ 200 lattice is the best option. 

To keep the updating procedure as simple as possible 
we used a standard Metropolis algorithm and the main 
runs where L = 100, (3 = 0.8,0.85 with 10 6 configura- 
tions, L = 200, (3 = 0.8,0.85,0.9 with 0.5 10 6 configura- 
tions. 

In Fig. 1 we present the effective exponent j\(y) for 
the L = 200 lattice at /3 = 0.80 and A = 1.0. As for 
other cases |6J the data points follow a smooth curve 
and a quadratic function fits perfectly the data. The 
extrapolated value for y — > is 7a(0) = 1.003(5): a sharp 
indication for a spontaneously broken CP symmetry at 
9 = tt. The same holds for different 0(1) values of A and 
for other lattice sizes and/or couplings. 

Following the procedure reported in ||| we can recon- 
struct, from the fits of y\{y), the order parameter de- 
pendence on 9. The results for the order parameter at 
= 0.80 and (3 = 0.85, for the L = 200 lattice, are shown 
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FIG. 2: Topological charge density versus 6/tt at /3 = 0.8, 
L — 200 (upper points) and j3 — 0.85, L — 200 (lower ones): 
first method (to guide the eyes the continuous curve joins 
the points) and from y\ extrapolation (diamonds, statistical 
errors are smaller than the symbol and are not reported). 



in Fig. 2. Indistinguishable results were obtained using 
the data of the L — 100 lattices, as expected if finite 
volume effects are negligible. 

On the same plot the output of the first method (based 
on the determination of the probability distribution func- 
tion) is also reported: the agreement between the two 
procedures is at (few) % level. These statements can be 
extended to the smallest gauge coupling available, with 
the only difference that vertical scale of the plot (and con- 
sequently the jump of the topological charge for 9 = tt) 
reduces increasing (3. For both methods the evaluation 
of statistical errors has been done computing the physi- 
cal quantities starting from independent simulation runs. 
It comes out that the statistical error on the topological 
charge density grows moderately with 9 and for 9 ~ tt 
takes its maximum value that, for the three values (3 = 
0.8, 0.85, 0.90 is of the order of 1, 2, 5% for the first 
method and 2, 3, 5% for the other one. 

The fitting procedure used in both methods are poten- 
tially sources of systematical errors. Concerning the first 
method |5| we remind that, in order to extract the prob- 
ability distribution function, we have to fit the h(x) data 
points. In principle any choice that respects the known 
symmetries (h(x) is an odd function) is admissible. We 
have tried different function families and the best results 
were obtained with h(x) = (a + bx 2 + ca; 4 )asinh(<ia;); 
with this 4-parameter function we were able to fit the 
data over the entire x interval we were interested in, 
with optimal y^/d.o.f. values. In particular we got 
X 2 /d.o.f. = 0.40 for (3 = 0.80, ^/d.o.j. = 0.62 for 
(3 = 0.85 and x 2 /d-o.f. = 0.67 for (3 = 0.90. 

For the second method 6J we also tried several different 
functions for fitting y\(y), obtaining fits of comparable 
quality, and identical results for the order parameter. 

Fig. 2 is a a posteriori check of the validity of our 
estimation of errors: since the two approaches are so dif- 
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FIG. 3: Topological charge density divided by the squared 
lattice size (in units of Al) versus 9/n for f3 = 0.80, 0.85, 0.90 
(continuous, dashed, dotted line respectively) in the L = 200 
lattice. 
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FIG. 4: As in Fig. 3 for the topological charge density divided 
by the measured topological susceptibility at 6 = 0. 

ferent each other, we do not expect agreement between 
the results if uncontrolled systematics plays a major role. 
This is, in our opinion, the most comfortable outcome of 
the whole procedure. 

The next question is to check scaling. In Fig. 3 we 
plot the continuum topological charge density x(9)/a 2 , 
in units of the lattice Al renormalization group pa- 
rameter, related to the lattice spacing a by the two 
loop perturbative renormalization group equation: a = 
A I ; 1 (27r/?) 2 / Ar exp(-27r/3). Note that a 2 varies by a factor 
of 4 in the considered (3 interval. In Fig. 4 we plot the 
lattice topological charge density normalized by the mea- 
sured topological susceptibility at 9 = 0. In both cases we 
plot the first method results only since, as stated above, 
those obtained with the second approach are practically 
not distinguishable. 

We see a good data collapse of plots at different (3 val- 



ues. In Fig. 3 violations to asymptotic scaling are at 
most ±6%. The origin of these violations could be the 
higher order corrections to the perturbative renormaliza- 
tion group equation and then, as seen in the plot, we 
would expect higher violations for the smaller (3 value. 
Indeed, as shown in Fig. 4, the non perturbative scal- 
ing is realized with much better accuracy: the three 
curves are now compatible within statistical errors. This 
is again a consistency check of both the results and the 
method(s). 

To conclude we want to emphasize that the contin- 
uum 6 dependence of a confining and asymptotically free 
quantum field theory has been fully reconstructed. Data 
collapse for different couplings within % level. The ev- 
idence for scaling in CP 9 model at non zero 9 is the 
strongest indication that the CP symmetry is sponta- 
neously broken in the continuum, as predicted by the 
large TV expansion. 

This work also shows how the inapplicability of the 
importance sampling technique to simulate physically in- 
teresting 9- vacuum models can be successfully overcome. 
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